# Wiring mache ich von links nach rechts !!!

...IDENTIFICATION DIVISION
...PROGRAM-ID FallingParticle
...VERSION 20240225
...COMMENT Calculate the deflection from the vertical caused by the Earth's rotation of a particle falling freely from rest from a height h.
...COMMENT Differential Equations:
...COMMENT    x''=-bz'+ay' # x-axis is along latitude,  directed to east
...COMMENT    y''=-ax'     # y-axis is along longitude, directed to north
...COMMENT    z''=-g+ax'   # z-axis is perpendicular to the surface of earth
...COMMENT    g: gravitational acceleration = 9,81 m/s²
...COMMENT    a: 2*omega*sin(phi)
...COMMENT    b: 2*omega*cos(phi)
...COMMENT    omega: rotation velocity of the earth = 2*pi/day
...COMMENT    phi: Latitude of location (0-90°)
...COMMENT    Initial Condition: z(0)=h
...COMMENT The full solution requires 6 INTEGRATORs, Anabrid-THAT just has 5. The deflection to longitude (y) is neglegible and can be omitted (marked #*).
...COMMENT It could also be solved in a separated algorithm omitting x.

...ENVIRONMENT DIVISION
...TIMEBASE 1ms

...DATA DIVISION
SET COEFFICIENT.1 TO AY # 2*2pi/day*sin(phi)
SET COEFFICIENT.2 TO AX # = AY
SET COEFFICIENT.3 TO B  # 2*2pi/day*cos(phi)
SET COEFFICIENT.4 TO G  # gravitational acceleration = 9,81 m/s²
SET COEFFICIENT.5 TO H  # height h
SET OUTPUT.X TO x
SET OUTPUT.Y TO y
SET OUTPUT.Z TO z

INITIALIZE H by -1 TO -h # same as COMPUTE -1 TIMES H TO -h
INITIALIZE G by +1 TO g

INTEGRATE 1*-bz',1*ay' TO -x' # Input is x''
INTEGRATE -x' TO x

INTEGRATE y'' TO -y'
INTEGRATE -y' TO y   

INTEGRATE z'' TO -z' 
INTEGRATE -z', IC:-h, LIMIT:(z >= 0) TO z

-x' * AX -> -ax'
-ax' = y''
-y' * AY -> -ay'
INVERT -ay' TO ay'
-z' * B -> -bz'
ADD 1*-ax', 1*g TO -g+ax'
-g+ax' = z''

...OPERATION DIVISION
...MODE REPEAT
...OP-TIME 7,3ms